Triple minima in Free Energy of Semifiexible Polymers 
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We study the free energy of the worm-like-chain model, in the constant-extension ensemble, as a function of the stiffness 
A for finite chains of length L. We find that the polymer properties obtained in this ensemble are qualitatively different from 
those obtained using constant- force ensembles. In particular we find that as we change the stiffness parameter, t = L/X, the 
polymer makes a transition from the flexible to the rigid phase and there is an intermediate regime of parameter values where 
the free energy has three minima and both phases are stable. This leads to interesting features in the force-extension curves. 
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The simplest model for describing semifiexible poly- 
mers without self-avoidance is the so called Worm-Like- 
Chain (WLC) model Jl]-^]. In this model the polymer is 
modeled as a continuous curve that can be specified by a 
d— dimensional (d > 1) vector x(s), s being the distance, 
measured along the length of the curve, from one fixed 
end. The energy of the WLC model is just the energy 
due to curvature and is given by 



(1) 



where u(s) — dx/ds is the tangent vector and satisfies 
u 2 = 1. The parameter k specifies the stiffness of the 
chain and is related to the persistence length A defined 
through (u(s) .u(s')) = e~\ s ~ s '/ A . It can be shown that 
K = {d- l)A/2 . 

The thermodynamic properties of such a chain can 
be obtained from the free energy which can be either 
the Helmholtz's (F) free energy or the Gibb's (G) en- 
ergy. In the former case one considers a polymer whose 
ends are kept at a fixed distance r [one end fixed at the 
origin and the other end at f — (0, ...0,r)] by an av- 
erage force (/) = dF(r, L)/dr, while in the latter case 
one fixes the force and the average extension is given by 
(r) = —dG(f, L)/df. It can be shown that in the thermo- 
dynamic limit L — > oo the two ensembles are equivalent 
and related by the usual Legendre transform G — F '— fr. 
For a system with finite L/X, the equivalence of the two 
ensembles is not guaranteed, especially when fluctuations 
become large. We note that real polymers come with 
a wide range of values of the parameter t = L/X [e.g. 
A w 0.1/im for DNA while A » lpm for Actin and their 
lengths can be varied] and fluctuations in r (or /) can be 
very large. Then the choice of the ensemble depends on 
the experimental conditions. Experiments on stretching 
polymers are usually performed by fixing one end of the 
polymer and attaching the other end to a bead which 
is then pulled by various means (magnetic, optical, me- 
chanical, etc.). In such experiments one can either fix the 
force on the bead and measure the average polymer ex- 



tension, or, one could constrain the bead's position and 
look at the average force on the polymer. In the former 
case, the Gibb's free energy is relevant while it is the 
Helmholtz in the second case. This point has been care- 
fully analyzed by Kreuzer and Payne in the context of 
atomic force microscope experiments Q. Theoretically, 
the constant-force ensemble is easier to treat, and infact 
an exact numerical solution has been obtained § (though 
only for t » 1). Data on force-extension experiments on 
DNA H have been explained using this ensemble || . The 
case of constant-extension ensemble turns out to be much 
harder and no exact solution is available. The t — ► and 
t — > oo cases correspond to the solvable limits of the 
hard rod and the Gaussian chain. The small and large t 
cases have been treated analytically by perturbation the- 
ory about these two limits @-|9|. Numerical simulations 
for different values of t have been reported by Wilhclm 
and Frey Jl(]], who have also obtained series expansions 
valid in the small t limit. A mean-field treatment has 
also recently been reported pi]] . 

In this letter we probe the nature of the transition from 
the Gaussian to the rigid rod with change of stiffness as 
shown by the form of the Helmholtz free energy of the 
WLC model (or equivalently the distribution of end-to- 
end distance). Extensive simulations are performed in 
two and three dimensions using the equivalence of the 
WLC model to a random walk with one-step memory. 
We find the surprising result that, over a range of values 
of t, the free energy has three minima. This is verified in 
a one-dimensional version of the model which is exactly 
solvable. 

We first note that the WLC model describes a parti- 
cle in d— dimensions moving with a constant speed (set 
to unity) and with a random acceleration. It is thus de- 
scribed by the propagator 



Z(x, u, L\x', u' ', 0) 



' (x' ,u') 
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where in the numerator only paths x{s), satisfying x(0) = 
x',x(L) — x, u(0) — u' and u(L) — u are considered. It 



1 



can be shown that the corresponding probability distri- 
bution W(x,u, L) satisfies the following Fokker- Planck 
equation |7|||- 







(3) 



where V? is the diffusion operator on the surface of the 
unit sphere in d— dimensions. The discretized version of 
this model is the freely rotating chain model (FRC) of 
semiflexible polymers 0. In the FRC one considers a 
polymer with N segments, each of length b = L/N. Suc- 
cessive segments are constrained to be at a fixed angle, 
8, with each other. The WLC model is obtained, in the 
limit 9, b -> 0, N -> oo keeping A = 2b/9 2 and L = Nb 
finite. 

Here we will consider the situation where the ends 
are kept at a fixed separation r [with x' at the origin 
and x — f — (0,...0, r)] but there is no constraint on 
u and u' and they are taken as uniformly distributed. 
Thus we will be interested in the distribution P(r, L) = 
(S(x — f)) = f duW(f, u, L): this gives the Helmholtz 
free energy F(r,L) = — Log[P(r, L)]. For the spherically 
symmetric situation we are considering, P(r, L) is sim- 
ply related to the radial probability distribution 5(7", L) 
through S(r,L) = Cr d ~ 1 P(r,L), C being a constant 
equal to the area of the d— dimensional unit sphere. It 
may be noted that the WLC Hamiltonian is equivalent 
to spin 0(d) models in one dimension in the limit of the 
exchange constant J — > oo (with Jb = n finite) and all re- 
sults can be translated into spin language. However, for 
spin systems, the present free energy is not very relevant 
since it corresponds to putting unnatural constraints on 
the magnetization vector. 

Numerical simulations: The simulations were per- 
formed by generating random configurations of the FRC 
model and computing the distribution of end-to-end dis- 
tances. To obtain equivalence with the WLC model the 
appropriate limits were taken. We note that because 
these simulations do not require equilibration, they are 
much faster than simulations on equivalent spin models 
and give better statistics. The number of configurations 
generated was around 10 8 for chains of size N = 10 3 . 
We verified that increasing N did not change the data 
significantly. As a check on our numerics we evaluated 



and 



Using Eq. (H) and following |12| we can 



compute these (in all dimensions): 
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Infact it is straightforward to compute all even mo- 
ments, though it becomes increasingly tedious to get 
the higher moments. Our numerics agrees with the ex- 



act results to around 0.1% for 



and 0.5% for 




FIG. 1. Monte-Carlo data for p(v,t) for the 

2— dimensional WLC for values of t = 10(o), 5, 3.33, 2 and 
1(V). The inset is a blowup of curves in the transition region 
(t — 4, 3.33, 2.86) and clearly shows the presence of the two 
maxima. Note that because of ±u symmetry, we have plotted 
data for positive v values only. For the fits at large and small 
t see text. 
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FIG. 2. Monte-Carlo data for p(v,t) for the 

3— dimensional WLC for values of t = 10(o), 5, 3.33, 2 and 
1(V). The inset is a blowup of curves in the transition region 
(t = 4, 3.85, 3.7) and shows the presence of the two maxima. 

The function P has the scaling form P(r, L) = 
jpp(r/L, L/X) and we will focus on determining the func- 
tion p(v,t) [|l|. In Fig. (|l|) and Fig. (||), we show the re- 
sults of our simulations in two and three dimensions. At 
large values of t there is a single maximum at?j = r/L = 
corresponding to a Gaussian distribution while at small 
t, the maximum is close to the fully extended value of 
v = ±1. The transition is first- order-like: as we decrease 
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t, at some critical value, p develops two additional max- 
ima at non-zero values of v. Further decreasing t weak- 
ens the maximum at v — until it finally disappears and 
there are just two maxima which correspond to the rigid 
chain. 

For the limiting cases of small and large values of t 
there are analytic results for the distribution function and 
as can be seen in Fig. (0J^) our data agrees with them. 
For large t we find that Daniels approximation |7j] , which 
is a perturbation about the Gaussian, fits the data quite 
well. In the other limit of small t the series solutions 
provided in |l0| fits our data. For intermediate values 
of t neither of the two forms are able to capture, even 
qualitatively, the features of the free energy. Specifically, 
we note that all the analytic theories(perturbative, series 
expansions and mean-field) predict a second-order-like 
transition and do not give triple minima of the free energy 
for any parameter value. 

It is instructive to study a one-dimensional version of 
the WLC which shows the same qualitative features (the 
equivalent spin problem is the Ising model). Consider a 
N step random walk, with step-size b which, with prob- 
ability e, reverses its direction of motion and with 1 — e, 
continues to move in the same direction. The appropri- 
ate scaling limit is: 6— >0,e^0,iV— >oo keeping 
L = Nb,t = LJ\ = 2Ne finite. Defining Z±(x, L) as the 
probability of the walker to be at x with either positive 
or negative velocity, we have the following Fokker-Planck 
equation: 



dZ± 
dL 



dx 



2A 



This can be solved for P(x, L) 
^p(x/L,L/X). We get 



= Z 4 



(5) 



Z_ = 



P (v,t) = — [ ±== +i (-VT^)\ 



(6) 



where Iq and I\ are modified Bessel functions. In Fig. || 
we have plotted p(v, t) for different values of stiffness. 
We find that it always has three peaks. Unlike in 2 and 
3 dimensions, the 5— function peaks at V = ±1 (which 
corresponds to fully extended chains) persist at all val- 
ues of stiffness though their weight decays exponentially. 
Similarly the peak at v = is always present. 

Discussion: The most interesting result of this pa- 
per is the triple minima in the Hclmholtz free energy 
of the WLC. Physically, this results from the compet- 
ing effects of entropy, which tries to pull in the poly- 
mer and the bending energy, which tries to extend it. 
This form of the free energy leads to a highly coun- 
terintuitive force- extension curve, very different from 
what one obtains from the constant force ensemble or 
from approximate theories. In Fig. (^) we show the 



force-extension curve for a two dimensional chain with 
t = 3.33. We see that there are two stable positions 
for which the force is zero. In the constant-force en- 
semble, it is easy to show that d(r)/df = (r 2 ) — (r) 2 
and so the force-extension is always monotonic. In the 
constant-extension ensemble, there is no analogous result 
(for finite systems) , and monotonicity is not guaranteed. 




FIG. 3. The exact distribution p(v,t) of the 

1— dimensional WLC [ Eq. (^)] for different values of t 
(10, 5, 3.33, 2, 1). Even for the most stiff chain considered here 
(t = 1), the distribution has a peak at the centre (in addition 
to the S— function peaks at ends) though it looks flat. 

Most of the recent experiments on stretching DNA 
have t > 100. The distribution is then sharply peaked 
at zero and one expects the equivalence of different en- 
sembles. Experimentally the value of t can be tuned 
by various means, for example, by changing the length 
of the polymer or the temperature. Polymer-stretching 
experiments can thus be performed for intermediate t 
values. Since we consider the tangent vectors at the 
polymer-ends to be unconstrained an accurate experi- 
mental realization of our set-up would be one in which 
both ends are attached to beads [see Fig. (pi)]. The 
beads are put in optical traps and so are free to rotate 
(this setup is identical to the one used in refn. E3|). 
Making the traps stiff corresponds to working in the 
constant-extension ensemble [3 and one can measure 
the average force. Our predictions can then be exper- 
imentally verified. We make some estimates on the ex- 
perimental requirements (for a 3-d polymer with stiff- 
ness t — 3.85). Assume that at one end, the origin, 
the trap is so stiff that the bead can only rotate. We 
make measurements at the other end. The trap-center is 
placed at fo = (0, 0, zq) and the mean bead displacement 
Az — {(z — zq)) gives the mean force (/) on the polymer. 
We then consider the problem of the polymer in the pres- 
ence of a trap potential V = [k t (x 2 + y 2 ) + k(z — z ) 2 ]/2. 
Assume k t >> k so we can neglect fluctuations in the 
transverse directions. The distribution of the bead's 



3 



position in the presence of the potential is given by 




FIG. 4. The free energy (dotted line) and the correspond- 
ing force-extension curve (solid line) for a 2— dimensional 
chain with t = 3.33. 
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FIG. 5. A schematic of the experimental set-up required 
to realize the constant-extension ensemble discussed in the 
paper (see refn. [14]). For a stiff trap the average displace- 
ment of the bead (Az) from the trap center is small and the 
average force on the polymer is: (/) = — k{Az). 

Q(f) = e-WW+nf)]/ J&fe-Pinti+vir)]. For a stiff 
trap, we can expand F about f = fo and find that the 
average displacement of the bead is given by: Az = 
Jd 3 f{z-z )Q(f) = -(f) /k', where k' = k + F"(z ) « k 
(valid except when zq ps L) and (/) = F'(zq). The rms 
fluctuation of the bead about the trap center is given by 
z 2 rms = k B T/k. Hence we get Az = -(f)z? m J(k B T) = 
-((f)L/k B T)(z 2 ms /L). The scaled force (f)L/(k B T) is 
of order 0.1. The different minima are separated by dis- 
tances 0.2L, hence to see the effect we need to have 
Zrms/L < 0.1. Thus finally we find that the typical dis- 
placement of the bead Az is about 0.01z rms . This is 
quite small and means that it is necessary to collect data 
on the bead position over long periods of time. 

As suggested in [|l0|, a more direct way of measuring 
the Helmholtz free energy would be to attach marker 
molecules at the ends of the polymer and determine the 



distribution of end-to-end distances. Fluorescence mi- 
croscopy as in [[l5| could be another possible method. 
It is to be remembered of course that real polymers are 
well-modeled by the WLC model provided we can neglect 
monomer-monomer interactions (steric, electrolytic etc.). 
Thus the experiments would really test the relevance of 
the WLC model in describing real semiflcxible polymers 
in different stiffness regimes. 

In conclusion we have presented some new and inter- 
esting properties of the WLC model and have pointed 
out that polymer properties are ensemble-dependent. In 
this paper we have given one example of qualitative dif- 
ferences in force-extension measurements in different en- 
sembles. Other quantitative differences will occur even 
in more flexible chains and should be easier to observe 
experimentally. We hope this work will motivate further 
experimental and theoretical work on this simplest model 
for semiflcxible polymers. 

We thank O. Narayan and J. Samuel for discussions. 
One of us (D.C.) thanks CSIR, India for support. 

Note: After submission of this paper, an exact numer- 
ical solution of the WLC model has been obtained, and 
has reproduced our results J16|. 
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